Surfactant in a Polyol–CO2 Mixture: Insights from a Classical Density Functional Theory Study

Silicone–polyether (SPE) surfactants, made of a polydimethyl-siloxane (PDMS) backbone and polyether branches, are commonly used as additives in the production of polymeric foams with improved properties. A key step in the production of polymeric foams is the nucleation of gas bubbles in the polymer matrix upon supersaturation of dissolved gas. However, the role of SPE surfactants in the nucleation of gas bubbles is not well understood. In this study, we use classical density functional theory to investigate the effect of an SPE surfactant on the nucleation of CO2 bubbles in a polyol foam formulation. We find that the addition of an SPE surfactant leads to a ∼3-fold decrease in the polyol–CO2 interfacial tension at the surfactant’s critical micelle concentration. Additionally, the surfactant is found to reduce the free energy barrier and affect the minimum free energy pathway (MFEP) associated with CO2 bubble nucleation. In the absence of a surfactant, a CO2-rich bubble nucleates from a homogeneous CO2-supersaturated polyol solution by following an MFEP characterized by a single nucleation barrier. Adding a surfactant results in a two-step nucleation process with reduced free energy barriers. The first barrier corresponds to the formation of a spherical aggregate with a liquid-like CO2 core. This spherical aggregate then grows into a CO2-rich bubble (spherical aggregate with a vapor-like CO2 core) of a critical size representing the second barrier. We hypothesize that the stronger affinity of CO2 for PDMS (than polyether) stabilizes the spherical aggregate with the liquid-like CO2 core, leading to a lower free energy barrier for CO2 bubble nucleation. Stabilization of such an aggregate during the early stages of the nucleation may lead to foams with more, smaller bubbles, which can improve their microstrustural features and insulating abilities.


INTRODUCTION
Polymer foams are lightweight materials with gaseous voids trapped in a polymer matrix. 1−7 Their properties depend strongly on microscopic features such as the size, density, and connectivity of the gaseous voids in the material. 8,9 The gaseous pores can be entirely separated from each other by the polymer matrix, leading to a foam with a closed-cell structure. 2,4,5,10 Alternatively, the gaseous pores can be interconnected within the polymer matrix, forming an opencell foam. 2,10 The closed-cell foams are rigid and are good thermal insulators. Consequently, they have found applications as materials in the construction, refrigeration, and automotive industries. Open-cell foams, on the contrary, are soft and flexible. These foams are used as materials for sound insulation and cushions for furniture, among many other applications.
Microscopic features of a foam are significantly influenced by the foam production process. 3 A standard procedure for producing foams 2,4 involves generating bubbles and stabilizing them within a polymer matrix. Reactive foaming takes advantage of chemical reactions between the blending reactants for gas evolution and their nucleation within a dense polymer medium. On the contrary, in a physical foaming process, the polymer is first saturated with gas at a desired pressure. Then the system conditions are instantly changed to initiate nucleation of gas bubbles in the system. This process yields a metastable condition in which the system is supersaturated with gas in the polymer. Such a supersaturated system evolves through nucleation of gas bubbles and subsequent growth of thus nucleated bubbles.
The presence of gas bubbles in a liquid makes a foam an inherently unstable system. 11,12 As a consequence, a foaming material ages over time. Drainage of liquid from the film between the bubbles, 13,14 bubble coarsening, 15 and bubble coalescence 16 are three main processes that contribute to foam instability. Surfactants are commonly used to stabilize foams against aging. 12,17,18 These molecules adsorb at the gas−liquid interface, improve interfacial properties, and constrain bubble coalescence and coarsening. 19−21 Many studies in the literature have investigated the effect of a surfactant on the insulating properties and the mechanical strength of surfactant-stabilized foams. 22−30 The most striking observation from these studies is that the addition of a surfactant yields a foam with a reduced cell size, an increased cell density, and an improved uniformity in cell size. Though the microstructural features of a foam are affected by both bubble nucleation and growth, Zhang et al. highlighted that they are more sensitive to the parameters governing bubble nucleation than those governing bubble growth. 31 While the role of a surfactant in stabilizing bubbles to achieve these properties has been reported, 32 its role in the nucleation of bubbles has not been reported, despite the high sensitivity of the nucleation barrier to interfacial tension. 33,34 We attempt to bridge that gap through this study.
Classical nucleaton theory (CNT), 35 classical density functional theory (cDFT), 36 and molecular simulation techniques 37 are commonly employed to investigate nucleation phenomena. CNT describes nucleation as formation of a new phase within a bulk phase. The formalism includes essential physics that governs the thermodynamic driving force for the formation of a new phase and the penalty for the formation of an interface between the new phase and the bulk phase. However, CNT uses equilibrium interfacial tension (IFT) and assumes a sharp interface. Such a theory is valid only near the coexistence. A semiempirical correction, called the Tolman correction, 38 is often employed to model the radius-dependent interfacial tension. When nanoscale bubbles are being studied, such a modification is less accurate. 39 There are also reports that invalidate the high free energy barriers predicted by CNT. 40 Additionally, the theory does not include any molecule specific aspects of the system. Molecular simulations can alleviate some of the problems posed by CNT. However, they suffer from finite size effects that arise while simulating nucleation in a small system with a fixed number of particles. 41 Though one can simulate a very large system 41 or simulate by imposing a constant chemical potential, 42 such simulations are computationally intensive.
cDFT 43,44 coupled with the string method 45 provides an appropriate framework for capturing the essential structure and thermodynamics associated with bubble nucleation. cDFT is a mean-field approach in which the free energy of the system is expressed as a function of spatially varying molecular densities. The density profile representing the equilibrium state of a system is determined by variational extremization of the free energy functional. If we have two such states, for example, state A being the homogeneous bulk and state B being a fully formed bubble of the prescribed size, a transition-state pathfinding algorithm like the string method 45 finds a minimum free energy path (MFEP) that connects states A and B. Though such an approach does not yield any dynamic information to go from state A to state B, it is still a very powerful technique for characterizing the MFEP and associated free energy barriers for the formation of a critical nucleus.
Xu et al. 46−49 have successfully employed cDFT with the string method to characterize the barriers for the nucleation of CO 2 bubbles in homopolymers like poly-methyl methacrylate (PMMA) and polystyrene (PS). In all of these studies, the authors have modeled their free energy functional for the cDFT based on the perturbed chain statistical associating fluid theory (PC-SAFT) equation of state (EoS). 50,51 Because the PC-SAFT EoS and cDFT have been demonstrated to quantitatively describe the gas solubility and interfacial properties in CO 2 −PS and CO 2 −PMMA systems, we employ the same approach to model surfactants and study their effect on bubble nucleation in polymer foams.
We specifically investigate the effect of a silicone−polyether (SPE) surfactant on CO 2 bubble nucleation in polyol. SPE surfactants are made of a polydimethylsiloxane (PDMS) backbone and polyether branches. 18,52,53 While alkylethoxylate surfactants 19 are ineffective in reducing the air−polyol interfacial tension, SPE surfactants are reported to significantly reduce the corresponding interfacial tension. 19 This makes SPE surfactants indispensable as stabilizers in the production of ubiquitous polyol-based foams like polyurethanes. The addition of co-solvents is known to promote or inhibit the surface-active abilities and aggregation behavior of surfactants in solution. 54−57 However, there have been no studies on how the presence of gas molecules would influence a surfactant's activity. Understanding these properties and how they manifest in the nucleation of gas bubbles is particularly relevant to characterizing the role of surfactants in foam production.
Using the cDFT based on the PC-SAFT EoS approach described above, we characterize the interfacial properties and aggregation behavior of SPE surfactants in a mixture of polyol and CO 2 . Then, using the string method, we compute the MFEP associated with the nucleation of a CO 2 bubble from a homogeneous mixture of CO 2 , polyol, and an SPE surfacant. Our main finding is that the SPE surfactant opens a low-energy barrier nucleation pathway. This has significant implications on the propensity for bubble nucleation and the resultant microstructural features of a foam.
The rest of the manuscript is organized as follows. We describe our models and cDFT approach in section 2. We report the results from our calculations and discuss them in section 3. We then conclude the article with an outlook on the path forward for developing foams with better physical properties.

Molecular Model and the Helmholtz Free Energy
Functional. We employ cDFT to model the CO 2 −polyol−surfactant ternary system. In our study, polyol is a homopolymer of poly(propylene oxide) (PPO) and the surfactant is a linear diblock copolymer with one block being PDMS and the other being PPO. Each of these molecule types is modeled as a tangentially connected chain of spherical beads. The Helmholtz free energy functional of such a system is then expressed as a sum of different perturbation contributions to the reference-state free energy functional.
To construct the cDFT, Xu et al. used weight-density functionals 58 and extended the PC-SAFT EoS to model the free energy functional for the inhomogeneous system. For the systems of interest to this work, we follow the same procedure and write the Helmholtz free energy functional [F({ρ})] as where F id is the ideal gas contribution, F hs is the excluded volume contribution due to hard sphere repulsion, F assoc is the association free energy due to the formation of a chain type molecule, and F disp represents dispersion interactions between segments of different molecule types in the system. Each of these perturbation contributions to the free energy density is expressed as a function of the spatially varying segmental density ({ρ}) of the different components in the system.
In the PC-SAFT EoS, the properties of a given molecule of type i are specified by the pure component parameters N i , σ i , and ϵ i , which represent the number of segments per molecule, the size of the segment, and the strength of the interaction between segments of the same type, respectively. The interaction potential between any two segments of types i and j is described by where σ ij = (σ i + σ j )/2 and = k (1 ) ij ij i j . k ij is the binary interaction correction term that is used to account for any missing interactions between segments of types i and j. If ρ i (r N i ) is the multidimensional density profile of molecule i with N i segments, where r N i = (r 1 , r 2 , ..., r ζ ..., r Nd i ), the corresponding segmental density [ρ i (r)] can be expressed as The different contributions to the Helmholtz free energy, listed in eq 1, are expressed as a function of these densities defined in eq 3. The ideal term (F id ) of the Helmholtz free energy is known exactly: where v is a volume scale that has no thermodynamic consequence as it just shifts the chemical potential by a constant. V B is the bonding potential that is used to enforce the chain connectivity between nearest-neighbor segments along the chain, and β = 1/k B T, where k B is the Boltzmann constant and T the temperature. V B is defined as The excluded volume contribution (F hs ), due to hard sphere repulsions, is modeled using the fundamental measure theory. 59 with where n j [{ρ}] = ∑ i n ji (j = 0, 1, 2, 3, V 1 , or V 2 ) are the Rosenfeld weighted density functionals. 58 These scalar and vector weighted density functionals are defined as In these equations, δ(r) is the Dirac delta function and Θ(r) is the Heaviside step function. The same weighted density functionals are used to describe other short-range interactions such as association and the local part of the dispersion interactions.
The thermodynamic perturbation theory of order 1 (TPT-1) 61−64 is used to model the contributions to the free energy due to the association type interactions (F assoc ). These interactions represent correlations within a molecule that arise due to the chain connectivity between the segments. If i = α, β, and γ index CO 2 , polyol, and surfactant segments, respectively, then the TPT-1 expression for F assoc is given by where N A and N B are the number of segments in the A and B type blocks of the surfactant chain, respectively, and N γ = N A + N B is the total number of segments per surfactant chain. Similarly, n 0,γ (r) = n 0,A (r) + n 0,B (r), where n 0,i (r) is the weighted density functional defined in eq 7. g ij (r) is the contact value of the correlation function between the segments of type i and j and is given by the expression in eq 9.
The contribution of dispersion (F disp ) to the free energy has local and nonlocal components. The local term (F disp−local ) is expressed as a perturbation to a chain-like reference fluid. 50,51,61,62,65,66 This is obtained by directly extending the corresponding PC-SAFT EoS expression to the inhomogeneous system using the weight density functionals defined in eq 7. The resulting expression is shown in eq 10.
where N̅ = ∑ i N i x i and x i is the mole fraction of molecule i. The constant coefficients, {a ln (m) |m = 1, 2; l = 0, 1, 2, ..., 6; n = 0, 1, 2}, are obtained by fitting the calculated binodal of the EoS with the experimental data for a great number of species. The values of these coefficients can be found in ref 51.
The F disp−local term in eq 10 alone is not sufficient to describe the contributions due to dispersion interactions. 67 A mean-field expression is included to account for any additional contributions Langmuir pubs.acs.org/Langmuir Article due to spatial inhomogeneity. 68 The corresponding nonlocal dispersion free energy term (F disp−nonlocal ) is given by

Euler−Lagrange Equations and the cDFT Numerical
Procedure. In a cDFT approach, the equilibrium state of the system is determined by minimizing its grand potential (Ω[{ρ}]). The Helmholtz free energy functional (F[{ρ}]) defined in eq 1 is related to the grand potential of the system by where μ i is the chemical potential of the ith molecule and ρ̂i(r N i ) and ρ i (r) are the corresponding molecular and segmental density profiles, respectively, as defined in eq 3.
, extremizing the grand potential with respect to the molecular densities [ρ̂i(r N i )] results in the following Euler−Lagrange equations.
Using eq 3 and the corresponding relation for the functional ], eq 13 can be re-expressed in Equation 14 can be further simplified by introducing a recursive function where the recursive function I ζ (r) is given by Equations 15 and 16 are the key equations for numerically computing the equilibrium density profiles of the different components in the system. In this work, we solve these cDFT equations in one-dimensional Cartesian coordinates and in spherical coordinates. The latter is used to study the properties of micellar aggregates and in the context of the string method to compute the MFEP for the nucleation of a spherical bubble.
The equations are solved in one-dimensional Cartesian coordinates to study the properties of the planar interface between the CO 2 -rich vapor and the polyol-rich liquid. To determine the spatial density profiles across the planar interface, we first compute the densities ( In our calculations, we choose L = 50σ 1 and Δz = 0.02σ 1 . Here σ 1 is the diameter of a polyol segment. Our choice for L is 8−10 times larger than the width of the interface. This allows us to accurately resolve the spatial density profiles at the interface and the bulk. The equations are solved using Picard iteration with a convergence criterion that the deviation of the Euclidean distance between two consecutive density profiles is <10 −6 . Because we solve cDFT equations in an open system, the location of the vapor−liquid interface is translationally invariant. Hence, we fix the spatial position of the interface. At the beginning of the iteration process, the position of the interface is at z = L/2 with the following spatial density profile: After every 10th iteration step, the entire profiles are shifted left or right, so that the Gibbs dividing surface 69 relative to the polyol is located at L/2.

Incipient Phase Calculation to Initiate Bubble Nucleation.
The primary objective of this work is to investigate the effect of the SPE surfactant on the nucleation of a spherical CO 2 bubble in polyol. To initiate bubble nucleation, we first saturate the polyol−surfactant mixture with CO 2 at the desired high pressure and temperature (303.8 K). We label this as the saturated state. The initial high pressure dictates CO 2 solubility in the saturated state. The desired CO 2 weight fraction in a foam formulation is ∼0.2−0.3 (w/ w), and this is realized at pressures of 6−7 MPa. For reference, the critical point in the CO 2 phase diagram is at 7.38 MPa and 303.8 K. Then, we instantly decrease the pressure to ambient conditions keeping the temperature and the CO 2 weight fraction in the system fixed. This leads to a metastable state in which CO 2 is supersaturated in the system. We refer to this state as the metastable parent phase.
We solve the PC-SAFT EoS at 1 atm pressure and 303.8 K to determine the densities of different components in the metastable parent phase, while keeping the CO 2 content unchanged from that in the saturated state. The metastable parent phase serves as a starting point for the nucleation of a CO 2 -rich bubble.
Because nucleation is a rare event, the system undergoes local density fluctuations representative of the microstates that are visited during the formation and breaking of subcritical nuclei. In our quasithermodynamic approach to nucleation, we seek to identify the CO 2rich bubble that is in chemical potential equilibrium with the metastable parent phase. We hypothesize that such a CO 2 -rich bubble is the incipient phase that the metastable parent phase tends to form; the composition of the incipient phase represents that of a large, wellformed bubble. We note that the pressure inside the incipient CO 2rich bubble is greater than the ambient pressure (i.e., the pressure of the metastable parent phase). As a result, the nucleated bubble eventually expands. In this work, we focus on the MFEP to the incipient CO 2 -rich bubble from a metastable parent phase and the surfactant effect on the associated free energy barriers. We solve cDFT equations with the string method 45,70 in spherical coordinates. At each point along the string, the system is at a constrained equilibrium, which allows the free energy of the bubble to be calculated. Connecting the points along the string results in the MFEP for the nucleation of an incipient CO 2 -rich spherical bubble.

String Method for the MFEP of Bubble Nucleation.
Within a mean-field framework, the MFEP for the nucleation is the most likely path that connects the initial and final metastable states via a transition state. 45 We solve eq 17 using a modified steepest descent algorithm to enforce the particular parametrization of the string. 48 The iteration starts with a set of initial density profiles between the initial state (s = 0) and the terminal state (s = 1, a well-developed bubble). States between s = 0 and s = 1 are obtained by linear interpolation. After each iteration, we reparametrize the states of the density profile equidistantly along the path. The process ends when the relative difference in the free energy along the path between two consecutive iterations is <10 −5 . As noted by Muller and co-workers, 71,72 such a steepest descent approach to constructing the string in the MFEP does not account for the dynamic constraint of local mass conservation. However, local extrema including saddle points in the free energy are unaffected by the local mass conservation. Therefore, we believe the qualitative picture for nucleation presented here remains valid.

Model Parameters and System Composition.
The pure component parameters for describing CO 2 and polyol come from the work of Xu et al. 47 and Ylitalo et al. 73 The corresponding parameters for describing different segments in the silicone polyether (SPE) surfactant are determined through a group contribution method. 74,75 These parameters are listed in Table 1. For a nonpolar or weakly polar system, such as the ternary system of interest to this work, the binary interaction correction term k ij that is used to correct for the dispersion interactions is negligibly small and is expected to have a minimal effect on the relevant properties of the system. Hence, in our modeling of ternary systems, we set the binary interaction correction terms to zero (i.e., k ij = 0).
For our study, we work with a model linear polyol with a molecular weight of 2700 g/mol. We choose to investigate the properties of different silicone surfactants whose overall molecular weights are close to that of the polyol. In this context, it is worth noting that a pure PDMS chain with a molecular weight as low as 1112 g/mol (equivalent to 40 PC-SAFT PDMS segments here) almost completely phase separates out of polyol. The presence of PPO type segments in the silicone surfactant is expected to improve the solubility of the silicone surfactant in polyol and render it surface-active. To systematically characterize the behavior of silicone surfactant in the CO 2 −polyol−silicone surfactant ternary system, we consider surfactants with different fractions of PDMS and PPO per chain. We acknowledge that the surfactant architecture affects their shape, aggregation, and interfacial properties. 76 However, in this work, we consider the simpler case of a PDMS Nf -PPO N(1−f) type linear diblock surfactant, where N is the total number of segments per chain and f is the fraction of PDMS segments per chain. SPE surfactants that are commonly used in foam formulations have f values of ∼0.3−0.6. 27 Hence, we restricted our analysis to surfactants with f values of ≤0.42. The different surfactants studied in this work are all listed in Table 2.
Such a choice of surfactants helps us to systematically investigate the trends in the aggregation behavior of a silicone surfactant in polyol, the effect of the CO 2 on it, and the vapor−liquid interfacial tension of the CO 2 −polyol−silicone surfactant ternary system. The results from these calculations are reported and discussed in the following section.

Micellization of a SPE Surfactant in Polyol.
The critical micelle concentration (CMC) is a characteristic property of a surfactant. It is the surfactant concentration in a solution above which most surfactant molecules exist in the form of aggregates. Micelles are reported to constrain the rate of drainage of the liquid from the film between bubbles, leading to a stable foam. 20 Micelles may also serve as seeds for gas adsorption, thereby promoting bubble formation through the heterogeneous nucleation pathway. Hence, the knowledge of the CMC of an SPE surfactant in polyol is crucial for the design of stable foams.
At any given surfactant concentration, micellar aggregates of different sizes form and break apart in the system. Their formation is always energetically favorable. However, when the surfactant concentration is lower than its CMC, the formation of these micelles (relative to the homogeneous bulk solution) is unfavorable due to the translational entropy loss of the individual surfactant molecules. If the surfactant concentration is higher than its CMC, micelles will form in large numbers. This serves as a working definition for determining the CMC of a surfactant. 77 In Figure 1a, we report the formation free energy of a micelle as a function of its size for different surfactant concentrations in the bulk solution. We define the micelle formation free energy (βΔΩ mf ) as β(Ω − Ω bulk ), where βΩ is the grand potential of the system containing the micelle and βΩ bulk is that of a bulk solution. A micellar aggregate of the desired size is obtained by restraining the polyol density, 78 to half its bulk value, at a given radial distance from the micelle center. To determine the size of such a micelle, we define the following order parameter: n exs−surf = ∫ dr[ρ surf (r) − ρ surf,blk ].   80 surfactant in polyol shown in Figure 1, the CMC is found to be ∼2.16% (w/w) surfactant in the polyol−surfactant mixture. For reference, the amount of surfactant that is commonly used in a foam formulation is ∼1−5% (w/w). A surfactant's tendency to aggregate into micelles can be enhanced by increasing the number of unfavorable interactions between the surfactant and the polyol. In this context, a key contributing parameter is the fraction ( f) of the PDMS segments of the surfactant. In Figure 1b, we report the micellization of different PDMS Nf -PPO N(1−f) type SPE surfactants in polyol. We fix N to 120 segments and vary f to study the effect of the length of the solvophobic block on surfactant micellization. We see that an increase in the fraction of PDMS segments from 0.25 to 0.42 leads to a decrease in the CMC from ∼15% (w/w) to ∼0.1% (w/w). With more solvophobic PDMS groups, the surfactant is less soluble in polyol. As a consequence, the surfactant tends to aggregate at lower concentrations.

Effect of CO 2 on the Micellization of the SPE Surfactant in
Polyol. In a foaming system, bubble nucleation is initiated by supersaturating the polyol−surfactant mixture with CO 2 . To investigate the effect of the surfactant on CO 2 bubble nucleation in polyol, it is then necessary to understand surfactant micellization as a function of the amount of dissolved CO 2 in the system. As noted in Models and Methods, we control the CO 2 content by the system pressure. The higher the system pressure, the more CO 2 is dissolved. In Figure 2, we report trends in the CMC of an SPE surfactant as a function of CO 2 content in the CO 2 −polyol−SPE surfactant solution.
When there is no CO 2 , the CMC of a model PDMS 40 -PPO 80 type SPE surfactant is found to be ∼2% (w/w) surfactant in the solution. An increase in the CO 2 content leads to a decrease in the CMC. When the CO 2 composition is ∼0.2% (w/w), the CMC is found to be as low as 0.1% (w/w) surfactant in the solution.
To understand the effect of CO 2 on the surfactant's CMC, we examine the spatial density profiles of PDMS and CO 2 segments in the most probable micellar aggregate (at the surfactant's CMC). We see from Figure 3a that the interface is located at roughly the same spatial position, irrespective of the CO 2 content of the solution. This suggests that all of these micelles have very similar micellar core sizes. However, an increase in the CO 2 solubility in the solution decreases the density of the PDMS segments in the micellar core. As shown     Figure 3b, this decrease in the density of PDMS segments in the micellar core is compensated by a significant increase in the CO 2 density in the micellar core. When a sufficiently large amount of CO 2 is dissolved in the mixture, the core of the micellar aggregate transforms from PDMS-dominated to CO 2 -dominated. Such a transformation is observed when the CO 2 concentration in the bulk solution is ≥0.2% (w/w). A cartoon representing the PDMS-dominated and CO 2 -dominated micellar aggregates is depicted in Figure  3c.
The preferential partitioning of CO 2 into the micellar core is a result of favorable CO 2 −PDMS interaction over CO 2 −polyol interaction. This results in swelling of the micellar core. However, we note that the micelle size is largely determined by the length of the PDMS block in the surfactant chain. As the micellar core is swollen due to the presence of CO 2 , the system now requires fewer surfactants to attain the same most probable micellar size. The presence of CO 2 in the micellar core may also relieve the stress due to the dense packing of the PDMS segments therein. This coupled effect due to the CO 2 may increase the drive for surfactants to aggregate more readily even at low surfactant concentrations.

Effect of the SPE Surfactant on the Vapor−Liquid Interfacial Tension in CO 2 −Polyol Mixtures.
A polymeric foam is characterized by the presence of an interface between the CO 2 -rich gas and dense polymer-rich medium. A surfactant's ability to reduce this interfacial tension facilitates the formation of such an interface. Here we investigate the interfacial tension between the CO 2 -rich vapor and polyol-rich liquid for the different surfactants studied in this work. The results of these calculations are summarized in Figure 4.
We find that, at a given weight percent of a PDMS f N -PPO (1−f)N type SPE surfactant in the solution, increasing the fraction of the PDMS segments in the surfactant leads to a steep decrease in the interfacial tension. This is a consequence of the lower CMC for the surfactants with a higher fraction of PDMS segments per chain. Intriguingly, the terminal interfacial tension, i.e., the interfacial tension at the CMC, is found to be nearly insensitive to the PDMS content per surfactant chain. This suggests that all of these surfactants have very similar surface-active abilities at their respective CMCs. However, for the same total surfactant weight percent in the solution, a surfactant with a longer block of PDMS is more effective at reducing the interfacial tension. It is also worth noting that these surfactants are moderate in their abilities to reduce interfacial tension, with an only 10−15 mN/m decrease in the interfacial tension before attaining the CMC. Very similar observations were recorded from the experimental investigations by Kendrick et al. 19 A surfactant's ability to reduce the interfacial tension between the CO 2 -rich vapor and the polyol-rich liquid informs us only of the reduced penalty to form such an interface. However, to make predictions about the rate at which bubbles are generated, we need to understand the effects of the surfactant on the bubble nucleation pathway and the associated free energy barriers.
In modeling SPE surfactant-mediated CO 2 bubble nucleation in polyol, we may consider nucleation when the surfactant's concentration in the solution is (a) below its CMC and (b) at or above its CMC. In the former scenario, one is likely to find polyol, CO 2 , and surfactant to be uniformly dispersed in the solution. When the system is brought into the metastable state by a sudden decrease in pressure, such a system is likely to undergo homogeneous nucleation. In the latter scenario though, the presence of preformed micellar aggregates presents complex heterogeneous pathways toward bubble generation. For example, (1) micelles may serve as seeds for bubbles to nucleate, (2) micelles themselves may evolve into bubbles, or (3) micelles may disintegrate into smaller aggregates and then evolve into bubbles. In principle, our models can be extended to investigate each of these pathways. As a first step toward understanding surfactantmediated CO 2 bubble nucleation, we investigate only scenario (a) in this study, i.e., homogeneous nucleation.

Effect of the SPE Surfactant on the Homogeneous Nucleation of CO 2 Bubbles in Polyol.
We show that the addition of the surfactant not only reduces the barrier for CO 2 bubble nucleation in polyol but also opens a new lower-energy barrier pathway through a spherical aggregate with a liquid-like CO 2 core.
In Figure 5a, we report the MFEP for CO 2 bubble nucleation in polyol when there is no surfactant. Here, the vertical axis represents the formation free energy of the bubble (relative to the homogeneous solution). The horizontal axis is an order parameter (V 2 ) that quantifies the size of the bubble. The definition for V 2 is We see that the MFEP between the homogeneous solution and a vapor-like CO 2 bubble passes through a single maximum in the free energy (see Figure 5a). At this maximum, the critical nucleus is CO 2 -rich and vapor-like, as seen in the density profiles in Figure 5c. Relative to the homogeneous solution (see Figure 5b), the free energy barrier corresponding to the critical nucleus is ∼35k B T.
Adding the surfactant leads to significant changes in the MFEP for bubble nucleation. The results of these studies are reported in Figure 6. First, the formation free energy of the critical nucleus decreases with the surfactant concentration.  Second, a shoulder appears during the early stages of the nucleation. This shoulder develops into a free energy barrier with an increase in the surfactant concentration.
To understand the importance of the free energy barrier that appears during the early stages of the nucleation, we choose the system with a ρ surf,bulk of 0.065% (w/w) as the model system and analyze the spatial density profiles along different stages of the nucleation. The relevant data are reported in Figure 7. Figure 7a shows that the first barrier has a formation free energy of ∼3k B T and is lower than the value of ∼11.8k B T that is noted for the second barrier. From the density profiles in Figure 7b, we note that the spherical aggregate representing the latter is characterized by a ρ COd 2 (r=0) of ≈0.6 g/mL. Due to the high density of CO 2 in the core, we label this aggregate as the one with the liquid-like CO 2 core. Formation of such an aggregate during the early stages of the nucleation and subsequent stabilization of the aggregate with a formation free energy of 0.61k B T (Figure 7c) may be a consequence of the favorable interactions of CO 2 −PDMS segments versus those of CO 2 −PPO segments. These aggregates eventually surpass an ∼11.0k B T barrier to vaporize into a spherical bubble with a vapor-like CO 2 core.
The lower nucleation energy barrier resulting from the opening of a two-stage nucleation pathway upon the addition of surfactant could yield higher nucleation rates compared to that of a surfactant free solution. Such a pathway could yield foams with more, smaller bubbles, which can improve their microstructural features and insulating abilities.

CONCLUSION
In this work, using cDFT, we investigated the effect of a silicone polyether surfactant on CO 2 bubble nucleation in polyol. We used a PC-SAFT EoS to model the chemically specific free energy functional for the cDFT calculations. Using these models, we first studied the interfacial and aggregation behavior of SPE surfactants in the CO 2 −polyol mixture. Following that, we computed the MFEP for CO 2 bubble nucleation in polyol and discussed the effect of an SPE surfactant on such an MFEP.
We find that the SPE surfactants aggregate into micelles in the CO 2 −polyol mixture. The terminal air−liquid interfacial tension of such a system is found to be 10−15 mN/m lower than that of the system without the surfactant. The CMC of an SPE surfactant is found to be strongly dependent on the CO 2 content of the system. The CMC decreases with increase in the CO 2 system. In a typical foam formulation, the surfactant's . CO 2 bubble nucleation in polyol when there is no surfactant in the system. (a) Minimum free energy path connecting the "homogeneous bulk solution" and a "spherical bubble" of predetermined size. V 2 is the order parameter directly related to the size of the spherical bubble. Spatial density profiles of different components in (b) the "homogeneous bulk system" (βΔΩ = 0) and (c) the system at the "critical nucleus" (βΔΩ = 35.23). Here, r is the radial distance from the center of the spherical bubble and βΔΩ = β(Ω s − Ω bulk ), where Ω s is the grand potential of the system with the bubble and Ω bulk is that of the homogeneous bulk solution. Figure 6. Minimum free energy path connecting the "homogeneous bulk solution" and a "spherical bubble of prescribed size" at different concentrations (ρ Surf,bulk ) of a PDMS 40 -PPO 80 type SPE surfactant in the system. Here, V 2 is the order parameter directly related to the size of the spherical bubble. For reference, the CMC of this SPE surfactant in CO 2 −polyol mixture is ∼0.1% (w/w). Langmuir pubs.acs.org/Langmuir Article composition is ∼1−5% (w/w) and CO 2 is saturated to ∼25% (w/w). Our calculations suggest that the surfactant's CMC is ∼0.1% (w/w) when the CO 2 content is ∼25% (w/w), implying that the surfactants may have aggregated into micelles even before the system is brought into the CO 2supersaturated metastable state. Similar to nanoparticles that impact foam production, 79,80 the presence of micelles may promote bubble nucleation, alter foam microstructure, and enhance foam stability. Though modeling heterogeneous nucleation is beyond the scope of this work, models similar to those presented in this work have been employed in the past to explore related problems like nanoparticle solvation in a polymer−CO 2 mixture. 81 In the future, we plan to extend the current model to study CO 2 bubble nucleation in the presence of preformed micelles. When the surfactant concentration in the CO 2 −polyol mixture is below its CMC, it is likely that the CO 2 bubble nucleates from a homogeneous solution of the CO 2 −polyol− surfactant mixture. In such a case, we find that the surfactant reduces the free energy barrier for CO 2 bubble nucleation in polyol. While the free energy barrier is ∼30k B T when there is no surfactant, it is found to be as low as 10k B T upon addition of 0.05% (w/w) SPE surfactant. Interestingly, the associated MFEP changes from a single-step nucleation process to a twostep nucleation process in the presence of a surfactant. From the density profiles, we find that the first barrier corresponds to the formation of a spherical aggregate with the liquid-like CO 2 core while the second barrier represents an aggregate with a vapor-like CO 2 core. We hypothesize that the formation of a spherical aggregate with a liquid-like CO 2 core during the early stages of bubble nucleation leads to a lower-energy barrier path for CO 2 bubble nucleation in polyol. This can enhance the nucleation rate and ultimately result in the production of foams with a reduced pore size and an increased number density of pores.
In this work, we focused on only the bubble nucleation aspects during the foam production process. However, for a comprehensive understanding, it is equally critical to investigate how surfactants constrain bubble coarsening and bubble coalescence. Stierle and Gross have recently reported a dynamic density functional theory (DDFT) to study bubble coalescence. 82 In their DDFT approach, the authors have accounted for the viscous forces as well as diffusive molecular transport through generalized Maxwell−Stefan diffusion. Similar to our work, Stierle and Gross use a PC-SAFT EoS to model the free energy functional for their DDFT. Suitable extension of their model to the SPE surfactant−CO 2 −polyol system would address the effects of surfactants on stabilizing polymer foams by constraining bubble coarsening and coalescence.